#include "opmatriz.h"
#include <stdio.h>

/*
 * inicia: aloca a matriz e retorna o ponteiro
 * */
PMatriz inicia(int nLin,int nCol){
	PMatriz a;
	(a = (PMatriz) malloc(sizeof(Matriz))) || die(MSG_ERR_MALLOC);
	a->nLin = nLin;
	a->nCol = nCol;
	int i,j;

	(a->elems = (double**) malloc(nLin*sizeof(double*))) || die(MSG_ERR_MALLOC);
	for (i=0; i<nLin; i++)
		(a->elems[i] = (double*) malloc(nCol*sizeof(double))) || die(MSG_ERR_MALLOC);

	return a;
}

/*
 * isVetor: retorna 1 se a matriz for um vetor
 */
int isVetor(PMatriz a){
	return (a!=NULL && (a->nLin == 1 || a->nCol == 1));
}

/*
 * isId: retorna 1 se a matriz é identidade
 */
int isId(PMatriz a){
	if (!isQuadrada(a)) return 0;

	int i,j;

	for (i=0; i<a->nLin; i++)
		if (!equals(a->elems[i][i],1.0)){
			return 0;
		}


	for (i=0; i<a->nLin-1; i++)
		for (j=i+1; j<a->nCol; j++)
			if (!equals(a->elems[i][j],0.0)){
				return 0;
			}


	for (i=1; i<a->nLin; i++)
			for (j=0; j<i; j++)
				if (!equals(a->elems[i][j],0.0)){
					return 0;
				}


	return 1;
}

/*
 * isQuadrada: retorna um se a matriz é quadrada
 */
int isQuadrada(PMatriz a){
	return (a!=NULL && (a->nLin == a->nCol));
}

/*
 * isEscalar: retorna 1 se a matriz é um escalar
 */
int isEscalar(PMatriz a){
	return (a!=NULL && (a->nLin == 1 && a->nCol == 1));
}

/*
 * getEscalar: se a matriz for escalar retorna esse valor
 * senão retorna menos infinito
 */
double getEscalar(PMatriz a){
	if (isEscalar(a)) return (a->elems[0][0]);
	else return (-HUGE_VAL);
}

/*
 * getColuna: (b) = coluna pi da matriz (a)
 * retorna 1 se conseguiu
 */
int getColuna(int pi, PMatriz a, PMatriz b){
	if (a==NULL||pi>=a->nCol) return 0;

	freeMatriz(b);
	int i;
	b = inicia(a->nLin,1);

	for (i=0; i<b->nLin; i++)
		b->elems[i][0] = a->elems[i][pi];

	return 1;
}

/*
 * leMatriz: le a matriz (a) nLin x nCol, do arquivo arq
 */
void leMatriz(FILE* arq, PMatriz a, int nLin, int nCol){
	freeMatriz(a);
	a = inicia(nLin,nCol);

	int i,j;
	for (i=0; i<a->nLin; i++)
		for (j=0; j<a->nCol; j++){
			if (arq==EOF) dieHelp(MSG_ERR_MF);
			fscanf(arq,"%lf ",&(a->elems[i][j]));
		}

}

/*
 * leVetor: le vetor (a) de dimensão dim
 * (coluna se isCol, linha caso contrário), do arquivo arq
 */
void leVetor(FILE* arq, PMatriz a, int dim, int isCol){
	freeMatriz(a);
	int nLin = isCol? dim:1;
	int nCol = isCol? 1:dim;
	a = inicia(nLin,nCol);

	int i,j;

	for (i=0; i<a->nLin; i++)
		for (j=0; j<a->nCol; j++){
			if (arq==EOF) dieHelp(MSG_ERR_MF);
			fscanf(arq,"%lf ",&(a->elems[i][j]));
		}
}

/*
 * clone: clona a matriz; (b) = (a)
 */
void clone(PMatriz a, PMatriz b){
	freeMatriz(b);
	if (a==NULL) return;
	int i,j;
	b = inicia(a->nLin,a->nCol);

	for (i=0; i<b->nLin; i++)
		for (j=0; j<b->nCol; j++)
			b->elems[i][j] = a->elems[i][j];
	}

/*
 * freeMatriz: desaloca a matriz (a)
 */
void freeMatriz(PMatriz a){
	if (a!=NULL){
	int i;
	for (i=0; i<a->nLin; i++){
		free(a->elems[i]);
		a->elems[i]=NULL;
	}
	free(a->elems);
	a->elems=NULL;
	free(a);
	a=NULL;
	}
}

/*
 * newId: (a) = matriz identidade dim x dim
 */
void newId(PMatriz a, int dim){
	freeMatriz(a);
	int i,j;
	a = inicia(dim,dim);

	for (i=0; i<dim; i++)
			a->elems[i][i] = 1;

		for (i=0; i<dim-1; i++)
			for (j=i+1; j<dim; j++)
				a->elems[i][j] = 0;

		for (i=1; i<dim; i++)
				for (j=0; j<i; j++)
					a->elems[i][j] = 0;


}

/*
 * impMatriz: imprime a martiz (a) no arquivo arq
 */
void impMatriz(FILE* arq, PMatriz a){
	if (a==NULL){ printf("pug\n"); return;}
	int i,j;
	for (i=0; i<a->nLin; i++){
		for (j=0; j<a->nCol; j++)
			fprintf(arq, "%f ",a->elems[i][j]);
		fprintf(arq, "\n");
	}

}

/*
 * impVetor: imprime vetor (a) no arquivo arq
 */
void impVetor(FILE* arq, PMatriz a){
	if (a==NULL) return;
		int i,j;
		for (i=0; i<a->nLin; i++)
			for (j=0; j<a->nCol; j++)
				fprintf(arq, "%f ",a->elems[i][j]);
		fprintf(arq,"\n");
}

/*
 * mult: multiplica duas matrizes; (c) = (a) X (b)
 * retorna 1 se consegui
 */
int mult(PMatriz a, PMatriz b, PMatriz c){
	if (a==NULL||b==NULL) return 0;
	if (a->nCol!=b->nLin) return 0;
	int i,j,k;

	freeMatriz(c);
	c = inicia(a->nLin,b->nCol);

	for (i=0; i<c->nLin; i++)
		for (j=0; j<c->nCol; j++){
			c->elems[i][j]=0;
			for (k=0; k<a->nCol; k++)
				c->elems[i][j] += (a->elems[i][k] * b->elems[k][j]);
		}
	return 1;
}

/*
 * transp: transpõe a matriz; (b) = (a)t
 */
void transp(PMatriz a, PMatriz b){
	freeMatriz(b);
	if (a==NULL) return;
	int i,j;
	b = inicia(a->nCol,b->nLin);

	for (i=0; i<b->nLin; i++)
		for (j=0; j<b->nCol; j++)
			b->elems[i][j] = a->elems[j][i];
}

/*
 * append: anexa duas matrizes lado a lado; (c) = (a)|(b)
 * retorna 1 se conseguiu
 */
int append(PMatriz a, PMatriz b, PMatriz c){
	if (a==NULL||b==NULL) return 0;
	if (a->nLin!=b->nLin) return 0;
	int i,j;
	freeMatriz(c);
	c = inicia(a->nLin,a->nCol+b->nCol);

	for (i=0;i<c->nLin;i++){
		for (j=0;j<a->nCol;j++)
			c->elems[i][j]=a->elems[i][j];
		for (j=0;j<b->nCol;j++)
			c->elems[i][j+(a->nCol)]=b->elems[i][j];
	}
	return 1;
}

/*
 * truncM: trunca a matriz (a) em nCol colunas, pondo em (b)
 * retorna 1 se conseguiu
 */
int truncM(PMatriz a, int nCol, PMatriz b){
	if (a==NULL||nCol<0||nCol>a->nCol) return 0;
		int i,j;
		freeMatriz(b);
		b =inicia(a->nLin,nCol);

		for (i=0;i<b->nLin;i++){
			for (j=0;j<b->nCol;j++)
				b->elems[i][j]=a->elems[i][j];
		}
		return 1;
}

/*
 * trocaCol: matriz (c) = (a) com a coluna colA trocada
 * pela coluna colB da matriz (b)
 * retorna 1 se conseguiu
 */
int trocaCol(PMatriz a, int colA, PMatriz b, int colB, PMatriz c){
	if (a==NULL||b==NULL) return 0;
	if (a->nLin!=b->nLin) return 0;
	if (colA<0||colA>=a->nCol) return 0;
	if (colB<0||colB>=b->nCol) return 0;
	clone(a,c);
	int i;
	for (i=0; i<c->nLin; i++)
		c->elems[i][colA] = b->elems[i][colB];
	return 1;
}

/*
 * inv: inverte a matriz (a), dado que só difere
 * uma coluna cCol de (b), e a inversa de (b) é conhecida
 * retorna 1 se conseguiu
 */
int inv(PMatriz a, PMatriz b, PMatriz invB, int cCol, PMatriz invA){
	if (a==NULL||b==NULL||invB==NULL) return 0;
	if (!isQuadrada(a)||!isQuadrada(b)||a->nLin!=b->nLin) return 0;
	if (cCol<0||cCol>=a->nCol) return 0;

	int i;
	PMatriz er;
	er = inicia(a->nLin,a->nCol);

	mult(invB,a,er);

	double vr = er->elems[cCol][cCol];
	if (equals(vr,0.0)) return 0;
	double fator = -(1/vr);

	PMatriz invEr;
	invEr = inicia(a->nLin,a->nCol);
	clone (er,invEr);

	for (i=0;i<cCol;i++)
		invEr->elems[i][cCol]*=fator;
	invEr->elems[cCol][cCol]=-fator;
	for (i=cCol+1;i<invEr->nLin;i++)
			invEr->elems[i][cCol]*=fator;


	PMatriz teste;
	teste = inicia(a->nLin,a->nCol);
	mult(invEr,er,teste);
	if(!isId(teste)) return 0;

	freeMatriz(invA);
	invA = inicia(a->nLin,a->nCol);
	mult(invEr,invB,invA);
	mult(invA,a,teste);
	if(!isId(teste)) return 0;

	freeMatriz(teste);
	freeMatriz(er);
	freeMatriz(invEr);

	return 1;

}




